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Abstract 

We consider classical TV-particle system with arbitrary central pair potential. Mechanical 
equilibrium condition in spherically-symmetric case leads to a nonlinear integro-differential 
equation for concentration n(r). For special state equation p = an 1 (1 + bn), original integro- 
differential equation transforms into purely integral one. Its solution (under 6 = 0) is written 
through a row over interaction parameter. Physical conditions for its convergency is discussed. 
For power-like potential kernal of integral operator is calculated in apparent kind. The cases 
of Coulomb and harmonic potentials are considered separately up to a third order. General 
scheme of application of the theory for some astrophysical and cosmological problems is pre- 
sented. Model system with spherical potential pits is considered. By perturbation theory first 
virial correction 8n is calculated. 

1 Introduction 

Imperfect systems, whose interparticle (quasiparticle) interactions principally determine its physi- 
cal properties, are of great interest today. Their distinctive peculiarities — nonlinearity, nonaddi- 
tivity of equilibrium thermodynamic potentials and parameters, critical phenomenas, selfordering 
effects etc - endow this sort of system by nontrivial dynamics and interesting observable effects. 
The last from the listed properties relates mainly to the imperfect systems with long-range inter- 
actions, which can provide orderization on scales about size of the system. Here and below the 
term " long-range" will mean order of magnitude of interactions radius equal or more then size of a 
system. The principal difficulty of statistical investigation of long-range interaction system (LRIS) 
concerns with a very complicated kind of statistical sum (integral) for long-range potentials. The 
well known virial decomposition jl| become invalid for LRIS, because of relation of interaction 
region to a volume of the system has order of unity. It was spontaneously found, that methods, 
relevant to the problem mainly are related to a different modifications of self-consistent field con- 
cepts. It is easily to show, that this is, in some sense, unavoidable. As it has been mentioned 
majority of customary properties of thermodynamic values are violated within LRIS. Particularly, 
well known rule for calculation of pressure: 

'-(Sf). 

becomes invalid or, more exactly, defines "generalized pressure". It is expressed by the sum: 
p = p^ +p( g \ where first term is ordinary local (contact) pressure, borned by classical or quantum 
motion of particles, while second one is determined by a long range interaction of some picked 
volume with environment. This decomposition corresponds to the subdividing of internal energy 
E = + E^ , where — local internal energy, — interaction energy of a small picked 

subsystem with environment. Now one can't understand p simply as a component of stress tensor 
of continuum media, if stress tensor of long-range field is not also included into consideration. The 
force components then will be 

Ft= { (°ik + 4k) d ^ k = f div + t ( % d ^ 
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where a® — ordinary continuum media stress tensor, t^' — field stress tensor, V — arbitrary 
volume inside the system. Mechanical equilibrium condition reads: 

div(trW+ttf>)=0, (2) 

that means intercompensation of local and global forces. But this resulting force depends on some 
internal properties of the system, for example — density of "generalized charge" distribution and 
we unavoidably go to the self-consistent problem. 

Attempts of investigations of LRIS have been made in a different areas: for dielectrics in 
electrostatic by thermodynamical methods in Q , for Newtonian gravity in balls stellar clusters by 
virial theorem JJ, by variational methods with fields equation as constraints — |4|, J5|, for quantum 
system — General problems are both choice of interparticle interaction or field equation, 

and solving of self-consistent equations, which are in general case nonlinear and integral. 

In present paper we don't restrict ourself to any concrete physical system and consider abstract 
centrally-interacting structureless classical particles. It is assumed, that pair potential satisfies 
some linear field equation, so that superposition principle is valid. We consider the case of identical 
particles, with pair potential: 

<f> = ag(r), (3) 

where a — generalized charge of a particle (the source of a long-range field ) , which is the same 
for all particles, r — distance between the particle-source and point of observation, g — universal 
for all particles function, depending on the distance and some universal constant parameters. It 
is postulated, that the system after a long time goes in equilibrium state, where all macroscopic 
parameters become constant. We are interested by practically important case of macroscopic 
spherically symmetry, when all macroscopic parameters (for example concentration) depends only 
on a distance from a center of the symmetry. 

Let in equilibrium state the system takes ball with radius R and let n(r) — concentration of 
particle, normalizing by the condition: 

R 

n dV = 4tt J r 2 n(r) dr = N, (4) 
v R o 

where N 3> 1 — full number of particles. In (^) central symmetry is used. 

Lets consider arbitrary point inside the system with radius- vector f. Potential dip at the point, 
producing by the particles inside some remoted volume dV in neighborhood of some another point 
r* of the system accordingly to (^) is equal: 

d<p(r)=g(\r-r*\)dL', (5) 

where dL' = a - n(r') • dV — total charge of dV'. Superposition principle for the full potential at 
r gives 

<p{r) = J g<y-?\)dL'. (6) 

Vr 



Taking into account apparent kind ^r— f*| = \/r 2 — 2rr' cos 9 + r' 2 and earring out angle integration, 
we get 



R 

I\I2 j I 



tp(r) = 2na j G(r,r')n(r')r' 2 dr' = 2naK[n] (7) 

where 



G(r, r') = / gi^/r 2 - 2rr'£ + r' 2 ) d£ = G(r' , r) (8) 



-l 

determines kernal Q of integral operator K. 

1 Rigorously speaking kernal in is G(r,r')r' 2 , but for the brevity we'll understand as kernal the function G. 
Note, that G(r',r) is, in fact, Green function for the linear equation, which pair potential ip obeys to. So by (^) we 
postulate special Green function, that corresponds to a wide class of linear differential field operators. 
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Attractive force, acting on a volume clement dV from the effective field of the system is 



dF = -^dL = • a • n(r) • dV, 
Or or 



and its volume density 



dF _ dp , . n , , , fdG(r,r') /2 



f = — = —zr-a ■ n(r) = -2na 2 n(r) / \-^-r' 2 n(r') dr' (9) 

dV dr w x ' J dr ' w 

o 

Equation (|^) of mechanical equilibrium, which in our case is reduced to —dp/dr + / = 0, with 
localQ pressure p, takes the form: 

R 

^ + 2na 2 n(r) [ dG( Z ,r 'K ' 2 n(r') dr' = (10) 
dr J or 

o 

The obtained equation belongs to a class of integro-differential ones. Really, pressure in general case 
is connected with concentration by thermal state equation p = p(n), which should be accounted 
by statistical methods. In our case it is determined also by long-range potential (its derivatives), 
whose influence on state equation is very complicated. 



2 Simplifying of equilibrium equation 

To make analytical investigation of eq. (|To|) some simplifying assumptions are necessary: 

1) Assume the local field has no influence on local state equation. This assumption will be valid 
in case of not strong fields and its gradients. In fact our consideration concerns with zero-th term 
of Taylor row of p over ip (or, more strictly its derivatives). 

2) Lets restrict our attention to a wide class of state equations of the form: 

p = mi 1 {l + bn) (11) 

where a, b and 7 — phenomenological constants. Note, that a has a sense of temperature multi- 
plier, which constantness within the system automatically implies its heat equilibrium. Expression 



(11) should be understood as the first two terms of effective virial expansion. For a not large 
concentrations second term much less then first one and so it influence on equilibrium distribution 
can be accounted by perturbations methods, with background, corresponding to the state equation 
p = an' . 

Stress, that thermal state equation ([ll]) is postulated by us, but not calculated as it should be 
done in fundamental theory. So our approach should be called semi- self- consistent: we take into 
account influence of long-range potential on concentration, but not on state equation. 

3) Lets take function g(r) from a ring of a generalized polynoms: 

<?(r)=5>r*S (12) 

i 

where e^, \ — arbitrary, but similar for all particles real parameters of potential. 
The first and second assumptions (under b = 0) in (|ll]) lead to equation: 

R 

07T1T-V + 2Tra 2 n(r) J ^^r' 2 n(r') dr' = 


or, after integration over r: 

n't- 1 = -sK[n] +A (13) 

in case 7 7^ 1, and 

n = Ae- sk[n \ (14) 

2 Here and below p will mean ordinary local pressure. 
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in case 7 = 1. Here A and A — yet undefined integration constants, s — 2ira 2 ('j — l)/a7, 

s = 2na 2 /a. (15) 

It easily to see, that equation ( |l4| ) is Boltzmann distribution written in self-consistent form. Really 
identifying (|l4| ) with n = noexp(—U / kT) under A = no, 

a = kT, (16) 

U = 2ira 2 K[n] (see (0), jl5|)), we obtain that this is true. 

The third simplifying assumption lets us by elementary integration calculate G. By linearity 
G over g it is sufficiently to carry out all calculations for one arbitrary power A in (|IJ). The result 
is: 

G(r, r>) = / e(r 2 2rr'£ + r' 2 f' 2 <% = [ , + f r + r> (17) 



A In 



3 General properties of equilibrium equations 



A 



In this section we present some general consequences of (|13j) and (|14[), which are independent of 
any assumptions concerning with perturbation theory (see following section). More exactly we are 
interested by possibility of cusps and finite radius R of the system, that is actual in galaxy and 
stellar dynamics fn], |§ |TJ, [TJ. The necessary conditions can be formulated by the two following 
simple theorems: 

Theorem 1. If g(r) is bounded from the up, then for 7 < 1 configuration with finite R is 
impossible. 

Theorem 2. // g(r) is bounded from the down, then for 7 > 1 configuration with cusp is 
impossible. 



Proof. Condition for finite R is n(R) = (R < 00). From (13) and ( |14| ) under 7 < 1 it is 
followed, that this possible if (and only if) K[n] — > +00. Bounded from the up potential g 
can be made negative on its whole domain of definition by finite shift: it only will redetermine 
integration constants A and A. Such negative g will give negative G and negative J^"[n], since 
n(r) > 0. So, .f£T[n], even it will have singular points will negative there: K[n] — > —00, that 

r—R 

means impossibility of finite radius. 

By the similar manner one can proof that under 7 > 1 and g bounded from the down K[n] at 

best gives: K[n] — * +00, that is opposite to that is required by cusps condition K[n] — ► 

— OO. 

The important consequence from the theorems is: finite configuration with cusp in case 7=1 
is possible, when g is unbounded both from the up, and from the down. This fact suggests, that the 
so called force softening procedure JT5[ , introducing in N-body simulation scheme to avoid local 
force infinities under collisions and to short computer time, can (under certain circumstances) 
drastically change global properties of investigated system. 

4 General form of solutions 

Let build at first a general solution for the case (p~4|) . For this purpose we write it through the 
formal row over powers of interaction parameter]^] s: 

n = X + X 1 s + X 2 s 2 + ..., (18) 

where Xi — yet unknown functions r. Expand then the exponent on K in ( |l4|) in the Taylors row: 

e-^M =1 _!^M + !!|!N + ... = (19) 



3 Inspite of dimension character of s we prefer it to dimensionless ordering parameter £ = a 2 R /a (see below), 
since R can be infinite. 
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s(K[X ]+SK[X 1 } + ...) , s 2 (K[X ] + sk[X 1 } + ...) 2 



1! 



2! 



To equate coefficients from the left and right in (|14|) after substituting (|18|) and (|19|), it is necessary 
to reexpand every term of ( p~Q| ) over s. The following useful formula, which can be proved by 
mathematical induction, is convenient: 



(x + xis + x 2 s 2 + . . .) p = Xq + 



(20) 



^s™ ^ p!(l!) fcl (2!) fe2 ■■■{n\) k 
Z^^T <L — = 



(k,n)—n 



(p - h - . . . - k n )\ 



k 1 k 2 ...k n X 



p—ki — k n k\ 



X, 



In this expression k = (hi, . . . , k n ), ft = (1, . . . , n) — integer-values n-dimensional vector, (a, b) 
- euclidian scalar product, coefficients B k k k can be calculated with the help of recurrent 
formulae, given in Appendix 1. For an integer positive p summation is carried out only for these 
k, which satisfies the nonequality p — k\ — . . . — k n > 0. It is easily to check that for particular 
cases p S j\f and p = — 1 under x. k — for all i = 2, . . . ( ppf ) reproduces well known Newton binom 
and formula for infinite geometrical progression correspondingly. 

Substituting ( |l8|) and ( [i"9| ) into (14), taking into account ( |20| ) and equating coefficients for 
equal powers s in (|14[) from the left and right, we get the system of recurrent relations on unknown 
coefficients in (|l8|): 

X = A; X 1 = -AK[Xo}; 



Xj = A 



3-1 



E 



(-i) p 



E 



(l!)^(2!) fc =---((i-p)!)^ , p 



(p- fci 



, (-l) j 



{k,j-p)=j-p 

(K[Xo]y- ki —- h i-*(£[Xi]) ki ■ ■ ■ c£pr,-_ p ])*'-> 

For the first four equations of the infinite chain ( |2l| ) gives: 
X = A; 

X x = AK[Xo\ = -A 2 K; 
X 2 = A{-K[Xi] + \K 2 [X \) = A 3 [K, 

x 3 = A(-k[x 2 ] + k[x ] ■ k[Xi] - ^k 3 [x })) = 

-A 4 [K 3 + \K[K 2 ] +K-K 2 + iif 3 ] 



(*[*<>] 



\K 2 ] 



(21) 



(22) 



Here K = K[l]; K n = K n [l]. 

The same way is appropriate for the more general case ([l3|) of equilibrium equation. Substi- 
tution of decomposition (|l^) (we omit a bar over letters) and equation of coefficients for all equal 
powers of s with (EOJ) gives: 



Xr'=A; 



JL_ v ( 7 -l)!(l!) fc i(2!)*»---(rc!) fc " 
n\ A (7-1-fci k n )\ 



(k,ft)—n 

or in resolved form 



7 — 1 — ^1 k n yk 



A -, • ■ ■ A r , 



X.,, 



(7-1) 



k[x n ^} 



- y 



( 7 -l)!(l!) fcl -((n-l)!)*-^ DT , 



(k,n—l)=n 



(7-1 fcn-l)! 



(23) 



7-1 fcn-1 i^fei 
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All noninteger factorials in ( p3| ) should be understood as infinite products a! = a(a — l)(a — 2) • 
The first four equation of the infinite chain has the form: 



A-r-i 



X, 



x 2 = 



7 — 1 



X, 



5-27 



(7-1)^ 
(7^1? 



(24) 



X 3 -^( 7 -2)^[K 2 



( 7 -2)(2 7 -3) g3 



(7 - 2)if • K 2 



At the end of the paragraph lets discuss some general physical considerations, relating to the 
convergency question of the row (jl^)[^ If generalized charges of a system will be too large or 
potential will be too rapidly fall in neighborhood of a particle (see following paragraph), then 
the configuration can not be equilibrium: there is singular condensate of the particles (with zero 
classical entropy) and row ( |l8[ ) is diverging. As a characteristic parameter, which is responsible 
for the such condensate forming, the value 



a 2 R> 



ur 
kT 1 



should be taken, where ur — potential energy of a pair particles at a distance of size of the system. 
The parameter £ is also characterize the influence of a long-range ordering forces, in comparison 
with influence of disordering heat chaotic particles motion inside the system. So, £ is natural to 
be called ordering parameter. 

On the other hand, an energy of long-range interaction of a particle with the whole system 
must be sufficiently large to provide bound state of the particles. In opposite case the system will 
be diffused. Corresponding parameters is 

a kT ' 



where Ur — potential energy of interaction of a particle with a whole system, 
physical conditions of equilibrium (and convergency of the rows) are: 



So reasonable 



£< KH. 



5 The cases A 



1 and A = 2 



Due to simplicity and practical importance particular interest have particles with Coulomb and 
harmonic potentials. Lets calculate the first nontrivial approximation X\ for 7 = 1. As it is seen 
from ( |22| ) the matter is reduced to a calculation of K. Simple analysis shows, that K — ► 00 under 
A < —3. In this case all terms of the formal row (|lj) become infinite, that expresses nonstability of 
configuration with these potentials: all particles falls onto the center. In case A > —3 elementary 
integration gives: 



K 



(A + 2)r U + 4 



A ^-2; 
(i? 2 



1 

27 



A = -2. 



1 



-((R + r) 



A+4 



(R 



A + 3 



-{{R + r) w + {R-r) A+i ) 



2Rr 



(25) 



4 Note, that direct mathematical investigation of the rows (pl|)-(|23[) convergency is problematic due to a very 
complicated structure of its coefficients. Even the evaluation of a number of nonzero B n for arbitrary fixed n, 
presented in Appendix 1, requires rather long calculations on n-dimensional integer-valued lattice. 
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For A = -1 we have K = eR 2 {l - z 2 /3). For A = 2, K = e(2i? 5 /5)(l + 5z 2 /3), 



r/R. In the 



For the former accordingly to (18) and (E2J) expression for n(r) takes the form: 



first case one should put e = —1, in second e = 1 to provide attractive character of interaction. 

(26) 



n ps riQ — UqsK = no + n sR (1 — 



which must vanish at z = 1. It gives uq = —3/2R 2 s. Substituting in (26), we get 

3 



4R 2 < 



1 



(27) 



Writing normalizing condition (Q) for (g^) we get the following first approximation for relation 
between equilibrium size of a system and number of particles: 



R C = ^N 
2tt 

Similarly for the case of harmonic pair potential we obtain: 

75 



128si? 5 



1-^ 



and 



Rh — 



5tt 
WsN' 



(28) 

(29) 
(30) 



Schematic curve for the both cases is shown in fig.|l](a). 

A n(r) (a) A L 



(b) 



r 



r 



R 



R 



Fi gUre 1: () — schematic dependency of n(r) for interparticlc Coulomb and harmonic potentials in the first approxima- 
tion; (b) — theoretical radial profile of luminosity for galaxies EO for the same potentials. 

The obtained expressions for equilibrium sizes of a system shows its strong nonadditivity, that 
has been mentioned in introduction. It is interesting, that equilibrium size of the system of bounded 
harmonic oscillators is decrease, when N increase, as it follows from (^p|). 

Note, that within central two-body problem there is interesting symmetry between this types of 
potentials: trajectories can be transformed to each other by the so called Bohlin's transformation 
fl6f| . As it can be seen from (^?j) and (|2^) the similar symmetry is revealed also within N-body 
problem at least in first approximation. As it follows from the theorems of the previous section, 
Coulomb system has R = oo and finite R in (^8|) is induced by approximation, while harmonic 
system has'nt cusp at the center. In fact, it means that Coulomb system of N < oo attractive 
particles totally diffuse in equilibrium. It theoretically support the fact, that galaxy and stars 
formations are looked as essentially nonequilibrium systems. Futher we'll use expressions (^7|)-(p9|) 
only to demonstrate the method, without any aims to reproduce observations. 

Using integrals K n from Appendix2, we can calculate the following approximation for n. The 
thiraQ has the form 



11 , 



1.26z 2 + 0.29^-0.03^ 



R : 



and for harmonic potential: 



0.52 

Irs 



l-0.95z 2 -0.14z 4 -0.05z 6 ) 



R 



sN 



0.24tt 
sN 



5 The second one is innormalizable on finite size R. 
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6 Application in cosmology 



Our model may be relevant to the systems of a large (N ^> 1) classical interacting objects, without 
any reference to its real scale. As an example, let consider cosmological application of the theory — 
EO galaxies or ball clusters. Statistical approach to such systems has been applied in [fn]. Rough 
consideration, based on virial theorem was made by Einstein in ||]. 

For a comparing the theory with observable datas it is necessary to translate it in terms of visual 
optical values, for example absolute surface luminosity. For a model, where all stars have the same 
overage characteristics it is easily to obtain the following formula, relating surface luminosity with 
concentration: 

\fR^f 

L = l J n(y/ ? + (31) 

where Iq — absolute luminosity of a single star, p — polar radius in the coordinate frame, related 
to the picture plane. 

Substituting of ( p7| ) or (^) in (|l]) leads to the following expression: 

L = A(l~?f\ C = |, (32) 

where 



coul " Its' aalm 32i?V 
which should be compared with experimental functions L(p). In figj^(b) the curve ( [32] ) is schemat- 
ically shown. This curve is well applied for a central part of cluster, where number of stars is 
sufficient to form self-consistent field, but individual stars are well resolved. Boundary region of a 
cluster is turn off from our consideration et al, since there are large fluctuations. The problem of 
separation and description of clusters core and bondary is considered in details in El . 

The method is also applicable for an investigation of a dark matter problem |17j |, which had 
been mentioned yet in Einstein cited work ||. The problem is inverse with respect to preceding 
one: using photometric datas one should reconstruct pair potential of interaction of astrophysical 
objects. Below we formulate general rough scheme of the approach, defering details and corrections 
for a future investigations. 

As it is well known from rotational curves analysis and from another astrophysical datas lumi- 
nous mass of a galaxy presents it smaller part, while remained one (about 90%) is confined in a 
dark halo, which extends to a distances about some radiuses of the visual central part. We assume, 
that equilibrium properties of a dark matter — its distribution and equilibrium size are practically 
independent from the luminous ones, so that the last can be viewed as probe subsystem at the dark 
background. Lets dark matter particles interact with each other by means of potential <p, where 
"tilde" will denote values, related to dark matter. Then using general formulas of our approach, 
one can calculate dark matter distribution h. This last create at every point additional potential 
for luminous matter, which is self interacting with Coulomb potential. The self-consistent "dark 
potential" ip, should be taken into account in expression (]?]) for a total potential: 



R 

12 



p(r) — 27ra J G(r, r')n(r')r' dr' + (p = 2naK[n] + (p. (33) 
o 

The potential (p when ft is known can be calculated with the help of (^), with interaction and 
temperature parameters a, a and kern G, taken in general form (|l^) with parameters A and e. 
These last are, in fact, to be found. For the above considered Coulomb and harmonic potentials 
simple calculations give: 



</?coul 



3a { z A z 2 1 \ . 15a /„ 7_ 9 \ 

TTT r + o H" COnst; <^harm = 7^7^ 1 + o Z + const > (34) 



4<5 V10 3 2/ ' r 2245 V 3 



where z = r / R, R — dark halo radius. Integration of equilibrium equation for luminous matter 
with potential (||^) in case 7 = 1 leads to the equation, generalizing (|l4|) (we omit bar over letters): 

n = Ae- sk ^- K *, (35) 



8 



where k = a' /a, a' — new interaction constant of luminous and dark matters. Introducing new 
function n* = ne K ^ , and new kernal G*(r,r') = G(r, r')e -Kl ^ r ), one can put (|35|) to the form (H): 



n* = Ae- sk ' [n " ] . (36) 

If luminous particles are interact with each other through Coulomb potential then with equality 
K* = e~ Ktp K , first approximation over s gives: 

« = -W)74:(l-4W) (37) 



AR 2 s V R 2 

— expression for luminous matter distribution with correction formfactors. These last are F±(ip) = 
(3 - 2e K *)~ 1 , F 2 = e-^/Fx = 3e~ 2KI ^ - 2e~ K *, where norming of potential is <p(R) = 0. Under 
K = (interaction between luminous and dark matter is absent) ( |37| ) transforms into (^7|). Then 
one should substitute expression ( pT^ ) into ( pl| ) for surface luminousity calculation. The result can 
be represented in the following convenient for comparison form: 

L(C) =L (C) + AL(C), 

where correction 



VH 5 

AL = ~H / (l + ^Wa-l)-^)^^^, (38) 



and Lo — surface luminosity (|32|) in absence of dark matter. For example, we present expression 
for AL, calculated on a dark matter background with Coulomb and harmonic potentials. Under 
the approximations R/R = s <C 1, which very simplifies final expression we have: 

ai _ A. - <T« (l - ic 2 ) ; 4- = ±; - ^. (39) 



where = (aa/aa). Comparison (|39| ) with observations gives us possibility for estimations of 
some combination of interaction parameters and kind of potential. Note, that proposed approach 
can can be complement to the so called MOND (Modified Newton Dynamics), where one deals 
with modified Newtonian pair potential too [pTs| . Obviously, the presented scheme is too rough 
and demands its following improving: account of a back reaction of a luminous matter on a dark 
one, more wide class of model potential, introducing of multicomponent systems etc.. 



7 The model with potential pits. 

Lets consider a system with potential of the following form: 

-1, r < p; 

(40) 

0, r > p, 

where p — length parameter, having the sense of effective interaction radius. The formulae (|1 
in this case is invalid, and direct calculation by (0) gives: 

^-(^^)- 

The integral K for the first approximation ( p^ ) has the form: 

1 ( R 2 , 2 R 2 , r 2 R 2 2rR 3 \ 
K =-2r{-2^-^)-— + —)- (42) 

Requirement of vanishing of expression n = hq — UqsK under r = R gives uq — — 4/sR(p 2 — R 2 /6). 
Under p 2 = R 2 /6 following approximations should be accounted. Expression for n after some 
algebra takes the form: 

r 2 R + r{p 2 -^R 2 )- R{p 2 - R 2 /2)\ . (43) 



sRr(p 2 - R 2 /6) 2 V 2 
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A n(r) (a) 




V2 



(b) 



> 



9/2 

Figure 2: (a) — schematic dependence of n(r) for the tips potential; (b) — dependency of equilibrium size of the system 
d = R/ p in effective interaction radius units on interaction parameter a — 3sN/8tt. 

Its schematic kind is shown in fig.|](a). Note, that in accordance with theorem 2, cusp — is induced 
by approximation. 

Nonnegativity condition for (|43|) will be nonequality: 



R < V2p, 



(44) 



that physically (in this approximation) means coincidence of orders magnitudes of size of the 
system and radius of pair interaction. Normalizing condition (Q) then gives: 



8^ R 2 p 2 
3s (p 2 ~R 2 /6) 2 



N. 



With ( |44| ) the (45) leads to the following restriction on a number N: 

sN < 12tt 



(45) 



(46) 



Physical sense of ( fig ) (which should remain valid, in principle, for any approximation) is concluded 
in the fact, that when interaction constant s is fixed, there is no equilibrium state for infinitely 
large N, since adding new particles we don't change attractive interaction intensity, while repulsive 
pressure is increased. This property simulate well known " saturation" phenomenas of water steam 
or nuclear systems. 

Dependency of equilibrium size d = Rj p on the interaction parameter a = 3sN/8ir, (0 < a < 
9/2) has the kind: 



d = 3V2^ 



1 



~\ 1+-CT 



and is shown in figj|(b). 



8 Virial correction calculation. 



Lets go backward to the expression ( |10[ ) and lets substitute in it the more general state equation 
(11). Using equality p' n = a , yn' 1 ~ 1 + 06(7 + l)n J , we get: 



dK 

ajn J ~ n + ab(-f + l)n 7 n' + 2wa n— — [n] = 

or 



(47) 



Dividing the both sides over n and integrating over r we go to equations generalizing ( |l3| ) and 
(0): 

-^—n'- 1 + afc^il?! 7 + 2Tra 2 K[n] + C = 0; (7/ 1); (48) 
7-1 7 

a\nn + 2abn + 2Tra 2 K + C = 0, (7 = 1). (49) 
We'll find solution to (^8|)-(p9|) by perturbation theory with small nonlinearity parameter b: 

n(r) = n(r) + b8\n + b 2 S2n . . . , (50) 
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where n — solution to G48|)-(f49|) under b — 0, which has been obtained earlier, 5\n,S-2n, 



unknown functions, defining corrections of different orders. Substituting ( |50| ) in (|4§|)-(49) and 
holding terms of a first order over b from the both side we get: 

dm + -^-n 2 + -K[Sin] = 0, 51 

— linear nonhomogeneous integral Fredholm equation on 5\n. This equation is valid for the case 
7 = 1 too. 

Lets find virial corrections for the considerd earlier potentials in case 7=1. 
a) Coulomb potential. Expression for the kernal G of K, as it follows from general formulae 
(17) has the form: 

Equation (^lj) with the kernal (|5|) and function n from (^7|) takes the form: 

+ i ( 1 - 4 V + ( 1 7TT ) {- ~ I r^n ^ ,lr' - 2 / r% n • ,lr \=il , ."W i 



8s 2 i? 4 V R 2 J AR 2 \ R 2 

It is convenient to go to the dimensionless variables z = r/R, y = 8\n ■ R e and to introduce 
dimensionless parameter £ 2 = R 2 /s 2 . Then equation (j5^) will take the following more simple 
form: 

V + le(l - z 2 ) 2 + |(1-^ 2 )|-| y z' 2 y ■dz'-ij z'y dz'\= 0. (54) 

Dividing over —3(1 — z 2 )/2z and differentiating twice over z, we get differential consequence of 
(54), which is differential equation of second order: 

u" + |(l-z 2 )«+|^ = 0, (55) 
where u — — 2zy/3(l — z 2 ). Its solution is given by the power row 

oo 

W (£,z) = ^> fc z fe , (56) 
where coefficients arc defined by the relations: 

a 2 = -^o, a 3 = - i (a 1 + 3£ ), a k+2 = 2(fc + 1)(fc + 2) , k = 2, . . . (57) 

Here ao and a\ — yet unknown integration constants. To find ones it is necessary to restrict 
obtained general solution to a subclass, satisfying (f54|). Substituting (56) in ( |54"| ) and integrating 
term by term, we equate coefficients for the same powers of z. This, apart from (|57]), will give 
additional equalities: 

3 °° a 

that fully defined our solution. Since ao = implies by ( Ejj) vanishing of all even coefficients in 
(56), we finally get: 

ftn = - l 2j?6 J E^ 2 "' where ( 58 ) 

*-!(^ + £<5T&ni>)- "> = >»+ A= 3( t(L" + 'i) 

From the ( |58| ) after some analysis we conclude that under b > 0, (that takes place for a 
majority of real systems and means additional repulsion) the correction bSin is always negative: 
system become more rared and its equilibrium size increased, as it is intuitively obviously. 
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b) Harmonic potential. Taking expression for n from ( |29|) and substituting it in ( |5l| ) we get 

y + 2fc 2 £ 2 (l - z 2 ) 2 + fc(l - z 2 )if [y] = 0, (59) 

where y and z — are the same, that in previous case, k — 75/128, £ = l/(si? 2 ). Formulae (|l7j ) for 
G in case A = 2 gives: 

G(r,r') = 2(r 2 +r' 2 ), (60) 

so for if [y] we get 



l 



= 2^ z' 2 (z 2 + z?)y(z!)dJ = Az 2 + B, 



where A — 2 J yz' 2 dz' , B = 2 J yz' A dz' . So the equation ( |59| ) on y is purely algebraic and its 



solution is reduced to finding of unknown A and B. To do it we substitute y from (J59|) in expressions 
for A and -B and after some transformations we go to the following linear system: 



(2fc/ 4 i + l)A + 2kI 2l B = -4k 2 $ 2 l 



22- 



2kI 61 A + [2khi + l)B = -Ak 2 i 2 h 2 ; 

where 



(61) 



/ ^-* 2 rd*=±(-ir k cz ' p+1 . (62) 



, : 2(g-fc)+P + l 

For the integrals in (|6l]) the (^2|) gives the following expressions: 

r _ 2 _ 8 2 8 2 

721 "15' 722 "105' 741 "35' 742 - 3 1 5 ' 761 "63- 
Substituting it into (^]J) and expressing A and B by Kramers rule, we get 

A = -J^Le » -0.024£ 2 ; B = _ J^L^ „ - .5lf; 

15808^ * ' 15808^ s ' 

Finally, 

0.39 

¥r 

As it has been in the previous case, the case b > rares the system. 



i ltl = — (l-2& 2 + l& 4 ). (63) 



9 Conclusion. 

So, the considered approach is relevant to the following problems: 

1) Calculating of equilibrium distribution of any LRIS; 

2) Finding of potential parameters by comparing calculated distribution with observation (in- 
verse to 1-st), that is more important from theoretical point of view, than 1-st. 

3) Calculation of virial corrections. 

The following obvious generalizations are possible: 

1) generalization on the case of coexisting chemical phases; 

2) generalization on rotating systems; 

3) generalization on identical quantum particles. This last case will involve nonlinear self- 
consistent multiparticle Shredinger equation. 

In conclusion I will be grateful to express many thanks to S.B.Moskowskiy for the interest and 
useful critical comments. 



A Appendix 1: coefficients 



Coefficients B% k , appearing in the rows (|20|) can be calculated for any n, hi, . . . k n with the help 
of the following recurrent formulae: 

1 ™ -2 l 

^fci...fc„_iO = ^ki—l k a ...fc„_i + £ (^J + l)^fci...fcj+lkj+i-l---fcn-i' 

J=1 (64) 

r>n i 

^OO.-.Ol ~~ 1 - 
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All another coefficients are zero, and beginning value B\ = 1. For a number f(n) of nonzero B n 
for a given n the following estimate can be obtained: 



f(n) < 



'2im 



Lets give the nonzero B^ k ^ for n = 1, 5 includingly. 
n = 1 = 1 



n = 2 B 2 



n — 3 Bfoo 



" = 4 B looo 



n — 5 -B| 0000 



i; 

= i; 
= 1; 



^2100 — 6; 



R3 
-°001 



R4 

-°1010 



i; 

= 4; 



^0200 ^ 3; 



R4 

-°0001 



B 5 

^31000 



B 5 1100 = 10; Bf 



00001 



1 ; -E>2o 100 
1 



1 o- b 5 

lu i -°12000 



15; B» 



10010 



5; 



(65) 



B Appendix 2: calculation of K n 



Lets show the results of calculations of integrals, which are necessary for third correction calcula- 
tion. 

For a Coulomb potential: 



K 2 = ^r- 



r 2 
3 


i? 2 ; 


r 4 


r 2 R 2 


30 


3 




R 2 r 


630 


30 



5i? 4 

~6~ ; 



+ 



5r 2 R 4 



61R 6 



18 



189 



i? 2 r 4 r 2 i? 4 



90 ' 
19i? 6 



27 



For harmonic potential: 



K x = 



Kn = 



2r 2 R 3 2R 5 



3 


+ 


5 ; 


8R 8 r 2 


+ 


184E 10 


15 


525 



704i? 13 r 2 768i? 15 

K 3 = + 



K[K 2 } = 



1575 2625 ' 

704i? 13 r 2 22336i? 15 



1575 



70875 
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